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Abstract. We present a mathematical description of the voter model dynamics on 
heterogeneous networks. When the average degree of the graph is [i < 2 the system 
reaches complete order exponentially fast. For \i > 2, a finite system falls, before it 
fully orders, in a quasistationary state in which the average density of active links (links 
between opposite-state nodes) in surviving runs is constant and equal to snj^b > while 
an infinite large system stays ad infinitum in a partially ordered stationary active state. 
The mean life time of the quasistationary state is proportional to the mean time to 
reach the fully ordered state T, which scales as T ~ fejjr^ ; where N is the number 
of nodes of the network, and fi2 is the second moment of the degree distribution. We 
find good agreement between these analytical results and numerical simulations on 
random networks with various degree distributions. 
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1. Introduction 

The voter model has become one of the most popular interacting particle systems [U [2] 
with applications to the study of diverse processes like opinion formation [3J H] , kinetics 
of heterogeneous catalysis [HIE], and species competition [Tj. The general version of the 
model considers a network formed by nodes holding either spin 1 or -1. In a single event, 
a randomly chosen node adopts the spin of one of its neighbors, also chosen at random. 
Beyond this standard version, several variations of the model have been considered 
in the literature, to account for zealots or inhomogeneities (individuals that favor one 
of the states) [8], constrained interactions [9], non-equivalent states [TO], asymmetric 
transitions or bias [11], noise [12] and ecological diversity [13]. It is also known that 
several models presenting a coarsening process without surface tension belong to the 
voter model universality class [14] . 

In a regular lattice, the mean magnetization, i.e., the normalized difference in the 
number of 1 and —1 spins, is conserved at each time step. Thus the magnetization 
is not a useful order parameter to study the ordering dynamics of the voter model. 
Instead, it is common in the physics literature to use as a order parameter the density 
of interfaces p, i.e, the fraction of links connecting neighbors with opposite spins. In a 
finite system, the only possible final state is the fully ordered state, in which all spins 
have the same value, either —1 or 1, and therefore all pair of neighbors are aligned 
(p = 0). These are absorbing configurations given that the system cannot escape from 
them once they are reached [15]. Despite its non-trivial dynamics, an exact solution 
has been obtained for regular lattices of general dimension d [6], becoming one of 
the few non-equilibrium models which are exactly solvable in any dimension. Indeed, 
the correspondence between the voter model and a system of coalescing random walkers 
helps to solve analytically many features of the dynamics [161 H7J- F° r d < 2, there is a 
coarsening process where the average size of ordered regions composed by sites holding 
the same spin continuously grows. In the thermodynamic limit, the approach to the 
final frozen configuration is characterized by the monotonic decrease in p, that decays 
as p ~ t~ x l 2 in Id and p ~ (hit) -1 in 2d [5]. For d > 2, the density of active interfaces 
behaves as p(t) ~ a — bt~ d ^ 2 [6J, thus p(t) reaches a constant value in the long time 
limit where the system reaches a stationary active state with nodes continuously flipping 
their spins. That is to say, full order is never reached. We need to clarify that the last 
is only true for infinite large systems, given that fluctuations in finite size lattices make 
the system to ultimately reach complete order. The level of order in the stationary state 
is quantified by the two-spin correlation function C^- = (SiSj) between spins i and j, 
that decays with their spatial separation r = \i — j\ as C(r) ~ r ( 2 ~ d ) [18], i.e, far apart 
spins become uncorrelated. Recent studies of the voter model on fractals with fractal 
dimension in the range (1,2), reveal that the system orders following p(t) ~ t~ a , with 
the exponent a in the range (0, 1) [T9J [20] . 

The voter model has recently been investigated on complex networks [2TJ [22j [23], [24J 
[25j [26j [27] , where its behavior seems to strongly depend on the topological characteristics 
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of the network. A peculiar aspect is that the dynamics can be slightly modified giving 
different dynamical scaling laws. For instance with node update, i.e., selecting first a 
node and then one of its neighbors, the conservation of the magnetization is not longer 
fulfilled. Instead the weighted magnetization is in this case conserved at each time step. 
With link update, where a link is selected at random and then one of its ends is updated 
according to the neighbor's spin, the conservation of the magnetization is restored [23J. 

A striking feature of the voter model on several complex networks, including Small- 
World, Barabasi- Albert, Erdos-Renyi, Exponential and Complete Graph is the lack of 
complete order in the thermodynamic limit. In this article, we provide an analytical 
insight of the incomplete ordering phenomenon in heterogeneous networks by studying 
the evolution and final state of the system using a simple mean-field approach. Despite 
that this approach is meant to work well in networks with arbitrary degree distributions 
but without node degree correlations, the qualitative results are rather general for 
many networks. We obtain analytical predictions for the density of active links (links 
connecting nodes with opposite spin) and the mean time to reach the ordered state as a 
function of the system size and the first and second moments of the degree distribution. 
These predictions explain numerical results reported in [231 [251 EHj and they agree with 
previous analytical results for ordering times [24|. 

The rest of the article is organized as follows. In section [21 we define the model 
and its updating rule on a general network. We then develop in section [3] a mean-field 
approach for the time evolution of the density active links and the link magnetization. 
This approximation reveals a transition at a critical value of the average connectivity 
fi = 2. When /i is smaller than 2, complete order is reached exponentially fast, whereas 
for fj, > 2, the system quickly settles in a quasistationary disordered state characterized 
by a constant density of active links whose value only depends on [i, independent on the 
degree distribution. We find that p is proportional to the product of the spin densities 
with a proportionality constant that depends on fi. This relation allows us to derive an 
approximate Fokker-Planck equation for the magnetization in section [U This equation 
is used in section [5] to study the relaxation of a finite system to the absorbing ordered 
state and in section [6] to obtain an expression for the survival probability of independent 
runs. The mean time to reach complete order, calculated in section U\ shows that the 
dependence of the results on the network topology enters through the first and the 
second moments of the degree distribution only. Convergence to the ordered state slows 
down as \i approaches 2, where ordering times seem to diverge faster than N. The 
summary and conclusions are provided in section [HJ In the appendix we present some 
details of calculations. 

2. The model 

We consider a network composed by a set of N nodes and the links connecting pair 
of nodes. We assume that the network has no degree correlations, i.e., the neighbors 
of each node are randomly selected from the entire set. We denote by Pk the degree 
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n active links k-n active links 




Figure 1. Update event in which a node i with spin Si = s (black circle) flips its spin 
to match its neighboring node spin Sj = —s (grey square). The possible values of the 
spins are s = ±1. Changes in the density of active links p and the link magnetization 
to = p ++ — p are denoted by Ap and Am respectively. 



distribution, which is the fraction of nodes with k links, subject to the normalization 
condition 'Yl lk Pk — 1- in the initial configuration, spins are assigned the values 1 or 
— 1 with probabilities given by the initial densities a + and cr_ respectively. In a single 
time step, a node i with spin Si and one of its neighbors j with spin Sj are chosen at 
random. Then i adopts j's spin (Si — > Si — Sj) (see Fig. [[]). This step is repeated until 
the system reaches complete order and it cannot longer evolve. 

3. Mean -Field theory 

In order to obtain an insight about the time evolution of the system we develop a mean- 
field (MF) approach. There are two types of links in the system, links between nodes 
with different spin or active links and links between nodes with the same spin or inert 
links. Given that a single spin-flip update happens only when an active link is chosen, 
it seems natural to consider the global density of active links p as a parameter that 
measures the level of activity in the system. 

In Fig. [H we describe the possible changes in p and their probabilities in a time step, 
when a node i with spin Si = s (s = 1 or —1) and degree k is chosen. We denote by n 
the number of active links connected to node i before the update. With probability n/k 
an active link (in this example i — j) is randomly chosen. Node i flips its state changing 
the state of its links from active to inert and vice- versa, and giving a local change of the 
number of active links An = k — 2n and a global density change Ap = 2 ■ Here 
pN/2 is the total number of links, p = (k) = ^2k ^Pk is the number of links per node 
or average degree. Assembling these factors, the change in the average density of active 



Analytical Solution of the Voter Model on Disordered Networks 



5 



links in a single time step of time interval dt = 1/JV is described by the master equation: 



dt ^ k dt 



, ^ 1 K ' ; k uN ' v ; 

fc ft fc ' n=0 ^ 

where -B(n, k) is the probability that n active links are connected to a node of degree k, 
and denotes the average change in p when a node of degree k is chosen. Given that, 
during the evolution, the densities of + and — spins are not the same, we expect that 
B(n, k) will depend on the spin of node i. For instance, when the system is about to 
reach the + fully ordered state, we expect a configuration where most of the neighbors 
of a given node (independent on its spin) have + spin, thus the probability that a link 
connected to a node with spin + (— ) is active will be close to zero (one). Therefore, we 
take B(n, k) as the average probability over the two types of spins 

B(n,k) = ya s B(n,k\s), (2) 

s=± 

where B(n, k\s) is the conditional probability that n of the k links connected to a node 
are active, given that the node has spin s. Replacing Eq. (j2J) into Eq. (TjQ) we obtain 

df> 2 yP k yo- s yB(n,k\s)l(k-2n) (3) 



(4) 



^ k s=± 



n=0 

(n) k ,s - ^(n 2 )k, s 



where (n)k, s and (n 2 )k, s are the first and the second moments of B(n, k\s) respectively. 

In order to calculate B(n,k\s) we assume that only correlations between the 
states of first neighbors are relevant, neglecting second or higher neighbors correlations. 
Therefore, we consider the conditional probability P{— s|s), that a neighbor of node 
i has spin — s given that i has spin s, to be independent of the other neighbors of i. 
This is known in the lattice models literature with the name of pair approximation, and 
it is supposed to work only in networks without degree correlations. Thus, B(n,k\s) 
becomes the binomial distribution with P(— s\s) as the single event probability that 
a link connected to % is active. P(— s\s) can be calculated as the average fraction of 
neighbors with spin — s to a node with spin s, i.e., the ratio between the total number 
ppN/2 of s —>■ — s links and the total number pa s N of links connected to nodes with 
spin s. We have used the symmetry in the states of the voter model and assumed that 
the average degrees of nodes holding spin 1 and —1 are the same and equal to /i. We 
have numerically checked that the last is valid for the original voter model, but if the 
two states are not equivalent or a biased is introduced, the average degrees are different. 
Then, P(—s\s) = p/2a s , and the first and the second moments of B(n, k\s) are 

kp 

(n) k ,s = ^ 

^ 1 * 2a, 4a? 
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Figure 2. Ratio between the density of active links and the product of the spin 

densities vs time in one realization of the voter model dynamics on networks with 

N = 10000 nodes and values of = 3,4,5,6,10 and 30 (bottom to top). Solid 

horizontal lines are the constant values 4 $ = 2 / A1 ~? . 

s (/i-i) 



Replacing these expressions for the moments in Eq. 
finally obtain 

dp 2p 



dt 



1 



P 



1 



and performing the sums we 



(5) 



2a + (l - <x + 

Equation ([5]) is the master equation for the time evolution of p as a function of the spin 
density o~ + (t). It has two stationary solutions, but depending on the value of /i, only 
one is stable. For p, < 2, the stable solution p = corresponds to a fully ordered frozen 
system. For /i > 2, the stable solution is 

p(t)=4Z{iM)a + (t)[l-a + (t)], (6) 
where we define 

(pi -2) 



(7) 



2(p-iy 

corresponding to a partially ordered system, composed by a fraction p > of active 
links, as long as <7+ ^ 0, 1. 

In Fig. [2] we test Eq. ([6]) by plotting the time evolution of the ratio between p 
and er + (l — a + ) in a single realization, for various values of \x. We observe that, even 
though the ratio varies over time, it fluctuates around the constant value 4£ predicted 
by Eq. ([6]). It is worth noting that the behavior of the ratio is the same from times 



Analytical Solution of the Voter Model on Disordered Networks 



7 



of order one to the end of the realization, where fluctuations increase in amplitude 
before the system reaches complete order. We also notice that fluctuations decrease as 
p increases, and they become zero in the complete graph case (p = N — 1), where we 
have pit) = 2a+(t) [1 - cr + (t)], for N > 1. 

In infinite large systems, fluctuations in cr+(t) vanish. Therefore, in a single 
realization we would see that cr + (t) = er+(0) for all t > and that p reaches an infinite 
long lived stationary state with p s = 4£cr + (0) [1 — cr + (0)]. Then, for networks with 
average degree p > 2, full order is never reached in the thermodynamic limit. 

In finite size networks, fluctuations eventually drive the system to one of the two 
absorbing states, cr + = 1 or cr + = 0, characterized by the absence of active links (p = 0). 
Although the parameter p is useful for finding an absorbing state, it does not allow 
us to know which of the two states is reached. For this reason we introduce the link 

magnetization m = p ++ — p , where p ++ (p ) are the density of links connecting two 

nodes with spins 1 (—1). It measures the level of order in the system, m = 1 (m = — 1) 
corresponding to the + (— ) fully ordered absorbing state and m = representing the 
totally mixed disordered state. Given that p becomes zero when m takes the values ±1, 
we guess that p should be proportional to 1 — m 2 . To prove this, we first relate a s with 
Pss by calculating the total number of links coming out from nodes with spin s. This 
number of links is pa s N, from which ppN/2 are s — > — s links, and p ss pN are s — > s 
links. We arrive to 

p ss = a s - p/2. 

Then, the link magnetization is simply the spin magnetization 

m = p ++ — p = o+ — er_ = 2a + — 1. (8) 

Combining Eqs. (jSJ) and (|HJ) we obtain that, neglecting fluctuations, p and m are related 
through the equation 

p(t) = e [1 - m\t)] . (9) 

Fig. [3] shows p vs m in one realization with p = 4 and iV = 10 4 . The system starts 
with equal density of + and — spins (m = and p = 1/2), and after an initial transient 
of order one, in which m stays close to zero and p decays to a value similar to £, p 
fluctuates around the parabola described by Eq. This particular trajectory ends at 
the (m = l,p = 0) absorbing state. 

4. Master Equation for the link magnetization 

In order to study the time evolution of the system we start by deriving a master equation 
for the probability P(m, t) that the system has link magnetization m at time t. In a time 
step, a node with spin s and degree k flips its spin with probability cr s P(— s\s) = p/2, 
after which the magnetization changes by Am = s5k, with 5^ = (see Fig. [1]), and 
with probability a s [1 — P{—s\s)] = a s (l — p/2a s ) its spin remains unchanged. We have 
used that the density of s spins and the conditional probability P{— s\s) in the subset 
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Figure 3. Trajectory of the system in a single realization plotted on the active links- 
link magnetization (p — rn) plane, for a network of size N = 10 4 and average degree 
p, = 4. Insets: Time evolution of m (left) and p (right) for the same realization. Wc 
note that p and m are not independent but fluctuate in coupled manner, following a 
parabolic trajectory described by p — -|(1 — m 2 ) from Eq.© (solid line). 



of nodes with degree k is independent on k and equal to the global density a s (this 
was first noticed in [21] and [26]). Using Eq. (Q we can write the probabilities of the 
possible changes in m due to the selection of a node of degree k as 

W m ^ m+5k = ^(l-m 2 )P k (10) 

W m ^ m = [1 - C (1 - m 2 )] P k . 

Thus, the problem is reduced to the motion of a symmetric random walk in the 
( — 1,1) interval, with absorbing boundaries at the ends and hopping distances and their 
probabilities that depend on the walker's position m and the degree distribution P k . 
The time evolution of P(m, t) is described by the master equation 

P(m,t + 6t)= J2 P *\ Wm 

P{m + 5 k , t) + W m _ Sk ^ m P(m- 5 k , t) 

k 

+ W m ^ m P(m,t)\ 
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k ^ 

1 [1 _ ( m _ 5 k f] P(m- S k , t) + [l- £(1 - m 2 )] P(m, t)|, 



where 5t = 1/N is the time step corresponding to a spin-flip attempt. In Eq. ffTTj) . 
the probability that the walker is at site m at time t + 5£ is written as the sum of the 
probabilities for all possible events that take the walker from a site m + A to site m, with 
A = 0, ±5 k and k > 0. The probability of a single event is the probability P{m + A, t) 
of being at site m + A at time t times the probability W m+/ \^ m of hopping to site m. 
Expanding Eq. ffTTj) to second order in m and first order in t we obtain 

d 2 P 



m 2 ) 



9m 



Thus, in the continuum limit (5t = 1/N — > as N — > oo), we arrive to the Fokker-Planck 
equation 

dP(m,t') d 2 r ,„ „. /s1 . 

where t' = t/r is a rescaled time, 

rs -^- = lj-'f (13) 

2 U/i)^2 (/i-2)/i 2 

is an intrinsic time scale of the system and /i 2 = Ylk k 2 Pk is the second moment of the 
degree distribution. We shall see in section [7] that the time to reach the ordered state 
equals r times a function of the initial magnetization. Note that, in complete graph, the 
corresponding Fokker-Planck equation derived for instance in has the same form 
as Eq. (TT2l) with t' = t/N, obtained as a particular case of a graph with distribution 
Pk — 3k,[i, n = N — I and /i 2 = /i 2 . The general solution to Eq. (fT2i) is given by the 
series expansion [281 129] 



P(m,t') = J2 A i C t /2 ( m ) e- {l+ ^ l+2 ^\ (14) 

1=0 

where A\ are coefficients determined by the initial condition and Cf^ 2 (x) are the 
Gegenbauer polynomials [20] page 980. Equation ffTTj) is of fundamental importance 
because it allows to find the two most relevant magnitudes in the voter model dynamics, 
namely, the average density of active links and the survival probability, as we shall see 
in sections |5] and [6] respectively. 



5. Approach to the final frozen state 

We are interested in how the average density of active links (p) decays to zero, where (•) 
denotes an average over many independent realizations of the dynamics starting from 
the same initial spin densities. Using Eq. (J9j) we can write 

(p(t')) =£(l-m 2 (t')> =ij dm{l-m 2 ) P{m,t'), (15) 
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with P(m,f) given by Eq. (1141) . The solution to the above integral with an initial 
magnetization m = 2<r + (0) — 1 is (see appendix Appendix A) 

(p(0> = £(1 - ml) e- 2t \ (16) 
and replacing back t' and £(//) we finally obtain 

(A* -2) , 



(p(*)> 



m je 



-2t/r 



17) 



2(/i-l) v 

We find that for p > 2, (p(t)) has an exponential decay with a time constant r/2, whose 
inverse gives the rate at which (p) decays. Given that r is proportional to A" (Eq. (Fl3|) ). 
the decay becomes slower for increasing system sizes. Eventually, in the limit of an 
infinite large network (p(t)) remains at the constant value £(1 — tuq) as it was discussed 
in section [3l while in a finite network, (p(t)) reaches zero in a time of order r. 

We have simulated the voter model on various types of random networks: degree- 
regular random graph (DR), Erdos-Renyi graph (ER), Exponential network (EN) and 
Barabasi- Albert network (BA). In Fig. H] we observe that the analytical prediction 
(Eq. (1171) ) is in good agreement with numerical simulations on these four networks. For 
a fix average degree p and system size N, r is determined by the second moment p2 of 
the network degree distribution P^. For these particular networks, p2 can be written 
as a function of p, because Pk only depends on p and k. As a consequence of this, 
r(p, N) is only a function of p and N. The values of r and p-i in the large A" limit are 
summarized in table [TJ For the case of DR, ER and EN, (p) is a function of t/N as it 
is shown in Fig. H] and pi is finite and independent on N. We have checked that the 
scaling works very well for networks of size A" > 100. For BA networks, pi diverges 



with A" (see calculation details in appendix Appendix B), thus we rescaled the x-axis 
by N/p2(N) in order to obtain an overlap for the curves of different system sizes. 



Network 


P k 


/J>2 


r(p,N) 


DR 






N 

(M-2) JV 


ER 




p{p+l) 


^(M- 1 ) AT 
( M +D(M-2) iV 


EN 




¥ 2 


4(^-1) AT 
5(M-2) J 


BA 




M(M+2) i_ / m(^+2) 3 N \ 


4 A1 (/.-l)Af/( / . 2 -4) 


2fc(fc+l)(fc+2) 




V (M+4) 4 1 


CG 


Sk,N-l 


(N — l) 2 


N 



Table 1. Node degree distribution Pk, its second moment \i2 and the decay time 
constant of the average density of active links r, for different networks. 



6. Survival probability 



In the last section we found that the density of active links, when averaged over many 
runs, decays exponentially fast to zero. In doing this average at a particular time t, 
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Figure 4. Time evolution of the average density of active links (p(t)) for (a) DR, 
(b) ER, (c) EN and (d) BA networks with average degree /i — 8. The open symbols 
correspond to networks of different sizes: N = 1000 (circles), N = 5000 (squares) 
and N = 10000 (diamonds). Solid lines arc the analytical predictions from Eq. (| 1 T|) . 
The average was taken over 1000 independent realizations, starting from a uniform 
distribution with magnetization mo = 0. 

we are considering all runs, even those that die before t and, therefore, contribute with 
p = to the average. In order to gain an insight about the evolution of a single run |25j, 
we consider the density of active links averaged only over surviving runs (p snrv {t)). If we 
define the survival probability S(t) as the probability that the system has not reached 
the fully ordered state up to time t, then we can write (p(t)) = S (t) (p snrv (t)) . 

In the Id random walk mapping that we discussed in section HJ S(t) corresponds to 
the probability that the RW is still alive at time t, that is to say, that it has not hit the 
absorbing boundaries m = ±1 up to time t. If at time t = 0, we launch many walkers 
from the same position m , each of which representing an individual run, then S(t) can 
be calculated as the fraction of surviving walkers at time t 

S(t) = j dmP(m,t). (18) 

The result of this integral for symmetric initial conditions (mo = 0) is given by the series 
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Figure 5. Survival probability S and average density of active links in surviving runs 
(p surv ) vs the rescaled time t/N for DR networks with degree ^ = 4 and sizes N = 100 
(circles), N — 400 (squares) and N — 1600 (diamonds). Top and bottom solid lines 
are the analytical solutions S(t) and (p surv ) = {pit)) / S(t) respectively obtained using 
equations (TIT))) and (TlT)) . 



(19) 



(see appendix |Appendix CD 

S{t) ~ h ^)n exp I ^ao ) 

As we observe in Fig.[5]there are two regimes. For t <C N, is S(t) ~ 1. For t > N/4, 
only the first term corresponding to the lowest I (I = 0) gives a significant contribution 
to the series, thus neglecting the terms with I > gives S(t) ~ |exp ^— T ^ N ^ j ■ For a 
general initial condition mo, we obtain that the survival probability decays as 
^3,, 2N / 2(u-2)u 2 * 



for t > N. (20) 



(A* - l)-" 2 N y 

Using Eqs. ( fTTj) and (!20l) we finally obtain that the density of active links in surviving 
runs is 

^-^(l-moV 2 ^ fort«AT; 

( f " % fort>N. 
3(/i-l) 
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Figure 6. Average height of the plateau for BA (circles) and DR (squares) networks 
of size N = 10000. The solid line is the analytical prediction -j^ryy- 



We find that the system reaches in a time of order N a partially ordered steady state, 
in which the average density of active links is 



(A* -2) 



(22) 



30* -1) 

In fig. [H] we plot the average height of the plateau as a function of \x obtained from 
numerical simulations on a Barabasi-Albert network and a degree-regular random graph. 
As Eq. ( 1221) shows, the average plateau value 2£/3 is only a function of the first moment 
of the distribution, as long as the network is random. The plateau is also independent 
on the initial condition mo, and the system size N for N large. 

A natural question is about the typical size of spin domains in the stationary state, 
where we use the term domain to identify a set of connected nodes with the same spin. 
Numerical simulations reveal that the system is always composed by two large domains 
with opposite spin until by fluctuations one of them takes over and the system freezes. 
This can be explained using percolation transition arguments on random graphs. Two 
connected nodes belong to the same domain if the link that connects them is inert, and 
this happens with probability q = 1 — p. Then, a domain that spans the system exists 
if q > q c = ^-j-, with k — — [21] ■ This gives a critical density 

2/i 



l>2 



\Ll-\L 



(23) 



Given that /i 2 > fJr, we have p c > ^— j- = 2£, and because the density of active links 
in one realization is equal or smaller than £ (see Fig. [3D, the system remains in the 
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"percolated phase", i.e., most of the nodes with the same spin are connected forming a 
giant domain of the order of the system size. 



7. Ordering time in finite systems 

A quantity of interest in the study of the voter model is the mean time to reach the 
fully ordered state when initially the system has magnetization m. In the random walk 
terminology of section HJ this is equivalent to the mean exit time T(m), i.e., the time 
that the walker takes to reach either absorbing boundary m = ±1 by the first time, 
starting from the position m. T(w) obeys the following recursion formula: 

T(m) = ^P fc [|(l-m 2 ) [T(rn + 6 k )+6t] 
k ^ 

+ |(i _ m 2 ) [T(m - 5 k ) + St] + [l - £(1 - m 2 )] [T{m) + St] j, 

with boundary conditions 

T(-l) = T(l) = 0. (24) 

The mean exit time starting from site m equals the probability of taking a step to a site 
m + A times the exit time starting from this site. We then have to sum over all possible 
steps A = 0, ±<5fc and add the time interval St of a single step. In the continuum limit 
(5k, St — > as iV — > oo) this equation becomes 

d 2 T(m) t 



dm 2 (1 — m 



2\ 



(25) 



where r is defined in Eq. (fT3j) . The solution to this equation is 
Tim) = t 



1 + m , / 1 + m \ 1 — m , f 1 — m 
In H In 



2 V 2 / 2 V 2 
or, in terms of the initial density of + spins a + = (1 + m)/2 

T(a + ) = - ^ - l ] f N [a + \na + + (l-a + )\n(l-a + )]. (26) 
(H - 2) fi2 

This expression differs with the one obtained in work [23] by a prefactor of 
However this factor does not seem to change the scaling of T(m) with the system size 
N, that was found to be in good agreement with numerical simulations. In Fig. (j7j) we 
show the ordering time t(o~ + ) as a function of the initial density of + spins, for a BA 
network with \i = 20, ER and DR networks with fi = 6. 

For a fixed N, Eq. (126!) predicts that T(m) diverges at \i = 2, but ordering times in 
the voter model are finite for finite sizes. To analyze this point, we numerically calculated 
T for an Erdos-Renyi network as function of \x for initial densities a + = o~_ = 1/2 (see 
Fig. [8]). For low values of \i, there is a fraction of nodes with zero degree that have 
no dynamics, thus we normalized T by the number of nodes N with degree larger 
than zero. As we observe in Fig. [SJ when /i decreases the analytical solution given 
by Eq. ( |26l) with fi 2 = p,(fi + 1) start to diverge from the numerical solution. This 
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Figure 7. Scaled ordering times vs initial density of + spins cr+ for networks of size 
N = 10 2 (circles), N = 10 3 (squares) and N = 10 4 (diamonds). Plots correspond to 
DR (top-left) and ER networks (top-right) with average degree fj, = 4 and BA networks 
(bottom- left) with /j, = 20. Solid lines are the analytical predictions from Eq. (j26|) . 

disagreement might be due to the fact that our mean-field approach assumes that the 
system is homogeneous, and neglects every sort of fluctuations, which are important in 
networks with low connectivity. However, we still find that T reaches a maximum at 
p ~ 2, where it seems to grow faster than N. 

8. Summary and conclusions 

In this article we have presented a mean-field approach over the density of active links 
that provides a description of the time evolution and final states of the voter model on 
heterogenous networks in both infinite and finite systems. The theory gives analytical 
results that are in good agreement with simulations of the model and also shows the 
connection between previous numerical and analytical results. The relation between the 
density of active links p and the density of + spins cr + expressed in Eq. (jHJ) allows to treat 
random graphs as complete graphs, and to find expressions for p and the mean ordering 
time in finite systems. For large average degree values, Eq. ([6]) reduces to the expression 
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Figure 8. Scaled ordering times vs average degree /i for Erdos-Renyi networks with 
N = 100 (circles), N — 1000 (squares) and N = 2000 (diamonds) nodes. The system 
size N was taken as the number of nodes in the network with degree larger than zero. 
The initial spin densities were c+ = cr_ = 1/2. The solid line is the solution given by 
Eq. 



for the density of active links in complete graph. Therefore, this work confirms that 
disordered networks with large enough connectivity are mean-field in character for the 
dynamics of the voter model. 

We find that when the average degree /i is smaller than 2, the system orders, while 
for jj, > 2, the average density of active links in surviving runs reaches a plateau of height 
3hi- % • ^ ue ^° fluctuations, a finite system always falls into an absorbing, fully-ordered 
state. The relaxation time T to the final absorbing state scales with the system size N 
and the first and second moments, \i and \ii respectively, of the degree distribution, as 

(A»-2)(U2 

Plateaus are also found on correlated networks with some level of node degree 
correlations, like for instance on small- world |2T, 26J, even though the plateau is lower 
than the one predicted by our theory. It might be interesting to modified the mean-field 
approach to account for degree correlations that correctly reproduce the behavior in 
very general networks. 
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Appendix A. Average density of active links 

To integrate Eq. (115j) . we use the series expansion Eq. (114ft for P(m,t') and write 

oo 

(p(0> = eE^ e ~ {m)( ' +2)t '> (a.i) 

1=0 

where we define the coefficient 

,2\ n^/ 2 



Di = j dm(l-m 2 )Ci /z (m). 

To obtain the coefficients Ai, we assume that the initial magnetization is m(t = 0) = m , 
i.e., P(m,t = 0) = 5(m — mo), from where the expansion for P(m,t') becomes 

oo 

A\ Cf^ 2 (m) = 5(m — m ). 

1=0 

Multiplying both sides of the above equation by (1 — m 2 ) Cf, (m) and integrating over 
m gives 

£ + ^ + ^ * *W = C 1 - -o) ^ /2 (-o) (A.2) 
i=o [ZL + d J 

where we used the orthogonality relation for the Gegenbauer polynomials Eq. MS 5.3.2 
(8) in page 983 of [30] with A = 3/2 



dmCf /2 (m) Cf/ 2 (m) (1 - m 2 ) = — — '-^-^—— 8^ (A.3) 



vrr(/ + 3) 



i 



All (Z + 3/2)[T(3/2)] 



and the identities T(l) = (I - 1)!, T(l + 1) = IT (I) and T(l/2) = v^F. Then, from 
Eq. ( 1A.2I) we obtain 

(2Z + 3)(l-m 2 )Cf /2 (m ) 

A ' " 2(1 + !)(/ + 2) • (A ' 4) 

To find D;, we use that the zeroth order polynomial is C^ 2 (m) = 1, together with 
the orthogonality relation Eq. ( 1A.3I) : 

D l = T dmCTim) cT(m) (1 - m 2 ) = —J^+^—S. 



-i 



4/! (/ + 3/2)[r(3/2)] 2W '° 



^2 1 +lK [ +2) 

(2/ + 3) ' V ; 

Then, using Eqns. (1A.4I) and (1A.5I) we find that the coefficients A\ and .D; are related 
by A\ Di — (1 — m 2 ,) Cf^ 2 (m ) ^,o- Replacing this relation in Eq. ( 1A.1I) and performing 
the summation we finally obtain 

<p(t')>=£(l-m 2 )e- 2i ', 

as quoted in Eq. ( |T6l) . 



Analytical Solution of the Voter Model on Disordered Networks 



19 



Appendix B. Calculation of fi 2 f° r Barabasi- Albert networks 

The Barabasi-Albert network is generated by starting with a number m of nodes, and 
adding, at each time step, a new node with m links that connect to m different nodes in 
the network. When the number of nodes in the system is N, the total number of links 
is mN, and therefore the average degree is [l — 2m. The expression for the resulting 
degree distribution, calculated for instance in [32], as a function of \i is 

P W = 2 t (ffl|(? + 2) - (B1) 



and its second moment is 

/i 2 = [ hmax ep(k)dk = ^ 1 "> I - — — — - (b.2) 



2 , //(// • 2) f Kmax kdk 



{k + l){k + 2) 



111 



2(k max + 2) 2 (n + 2 



_ (k max + l)(fi + A) 2 _ 

The lower limit fi/2 of the above integrals correspond to the lowest possible degree m, 
since nodes already have m links when they are added to the network. The reason for 
an upper limit k max is that the contribution to /i 2 from large degree terms is important 
due to the slow asymptotic decay P(k) ~ k~ 3 , unlike for instance in Erdos-Renyi or 
Exponential networks where P(k) decays faster than fc~ 3 , thus high degree terms become 
irrelevant. k max is estimated as the degree for which the number of nodes with degree 
larger than k max is less than one. Then 

1 + 2) f°° dk _M/i + 2) ln / (k max + l) 2 



max ' 



N 2 J kmax k(k + l)(k + 2) 4 \k max (k max + 2) i 

Assuming k max 3> 1, the expansion of the logarithm to first order in l/k max is 1/k 2 
Then, solving for k max , we obtain 

k max ~ v /u( M + 2)/4iV 1 / 2 , (B.3) 

i.e, the maximum degree diverges with the system size. 

Taking k max ^> 1 in Eq. ( 1B.2[) and replacing the value of k max from Eq. ( IB. 31) gives 
the expression quoted in table [H for the second moment of a BA network 

, 2= ^, n ^ + 2 )^V (B.4) 



(f + i) 

Appendix C. Survival probability 

By using the series representation Eq. (JUJ), the survival probability quoted in Eq. (jT8l) 
can be written as 

oo 

S{t) = Y,AiBie^ l+1){l+2)t ' ) (C.l) 

1=0 

where we define 



B t = dmCf /2 {m). (C.2) 
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To obtain the coefficients B[, we use the derivative identity Cf^ 2 (m) = 
derived from Eq. MS 5.3.2 (1) in page 983 of [30J with A = 3/2. Then 

B, = CtlKl) - O-l) = 1 - (-!)« = {° \ ^ (C.3) 

where we have used the relations C^ 2 {1) = 1 V/ and C^ 2 {— 1) = (—1)' that follow 
from Eq. MO 98 (4) (page 983) and the parity of the polynomials (page 980) of [3D] 
respectively. 

An explicit function for the coefficients A\ of Eq. ()A.4j) can only be found for the 
mo = case, given that for mo 7^ it seems that a closed expression for the polynomials 
Cf^ 2 (mo) cannot be obtained. To obtain the coefficients Cf^ 2 (0) we use the recursion 
relation Eq. Mo 98 (4) (page 981) of for m = x = and A = 3/2, together with the 
values of the zeroth and first order polynomials Cg^ 2 (0) = 1 and Cf^ 2 (0) = 0. Then 

r 3/2 fm _ ([+!) r 3/2 rn x _ / /Odd 

C, (0) - — 6 ^ 2 (°) - I ( _ iy /2 z even (C-4) 

Plugging the above expression into Eq. (1A.4D gives A; = for / odd and 

Ai — 2(;+2)ii even. 

Then, using Eq. (IC.3I) . the product A t Bi can be written as 

[O /odd 
AiB t =l { -iy/2 {2 i + 3)(i-i)v. j (C5) 

^ (J+2)!! ' e e 

Finally, making the variable change / — > 2 /, Eq. (IC.lj) becomes 

^=g: ( - 1) '%t h 1 ( ) r i) " e " 2 ' 2,+1)( ' +1) '' (c - 6) 



Z=0 



Replacing t' by £/t(/x, iV), we obtain the expression quoted in Eq. (fl9l) . 



